No Haiti data

rmshare.sum.6.noha<-rmshare.sum.6[!rownames(rmshare.sum.6) %in% "ha",]
rm.nebf<-metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=rmshare.sum.6.noha,sm="RD", backtransf=FALSE)
forest(rm.nebf,smlab="Standardized \n microbiome age difference",sortvar=rmshare.sum.6.noha$pop)

rm.nebf

                                          RD            95%-CI %W(fixed)
Thompson et al 2015 ( USA(NC) )       0.8754 [ 0.1366; 1.6142]       1.5
Wood et al 2017 ( South Africa )      0.1214 [-0.2963; 0.5390]       4.8
Pannaraj et al 2017 ( USA(CA_FL) )    0.2665 [ 0.0464; 0.4866]      17.4
Subramanian et al 2014 ( Bangladesh ) 0.0568 [-0.0605; 0.1741]      61.4
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.7355 [ 0.3600; 1.1109]       6.0
Azad et al 2015 ( Canada )            0.6149 [ 0.3048; 0.9250]       8.8
                                      %W(random)
Thompson et al 2015 ( USA(NC) )              8.1
Wood et al 2017 ( South Africa )            14.8
Pannaraj et al 2017 ( USA(CA_FL) )          20.4
Subramanian et al 2014 ( Bangladesh )       22.9
Sordillo et al 2017 ( USA(CA_MA_MO) )       15.9
Azad et al 2015 ( Canada )                  17.8

Number of studies combined: k = 6

                         RD           95%-CI    z  p-value
Fixed effect model   0.1988 [0.1069; 0.2908] 4.24 < 0.0001
Random effects model 0.3835 [0.1247; 0.6422] 2.90   0.0037

Quantifying heterogeneity:
tau^2 = 0.0726; H = 2.20 [1.49; 3.25]; I^2 = 79.3% [54.7%; 90.5%]

Test of heterogeneity:
     Q d.f. p-value
 24.11    5  0.0002

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=rm.nebf$studlab,pval=rm.nebf$pval))

No USA(NC) data

rmshare.sum.6.nounc<-rmshare.sum.6[!rownames(rmshare.sum.6) %in% "unc",]
rm.nebf<-metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=rmshare.sum.6.nounc,sm="RD", backtransf=FALSE)
forest(rm.nebf,smlab="Standardized \n microbiome age difference",sortvar=rmshare.sum.6.nounc$pop)

rm.nebf

                                           RD            95%-CI %W(fixed)
Wood et al 2017 ( South Africa )       0.1214 [-0.2963; 0.5390]       4.8
Pannaraj et al 2017 ( USA(CA_FL) )     0.2665 [ 0.0464; 0.4866]      17.4
Bender et al 2016 ( Haiti )           -0.1769 [-0.8203; 0.4664]       2.0
Subramanian et al 2014 ( Bangladesh )  0.0568 [-0.0605; 0.1741]      61.1
Sordillo et al 2017 ( USA(CA_MA_MO) )  0.7355 [ 0.3600; 1.1109]       6.0
Azad et al 2015 ( Canada )             0.6149 [ 0.3048; 0.9250]       8.7
                                      %W(random)
Wood et al 2017 ( South Africa )            14.3
Pannaraj et al 2017 ( USA(CA_FL) )          20.4
Bender et al 2016 ( Haiti )                  9.1
Subramanian et al 2014 ( Bangladesh )       23.1
Sordillo et al 2017 ( USA(CA_MA_MO) )       15.5
Azad et al 2015 ( Canada )                  17.6

Number of studies combined: k = 6

                         RD           95%-CI    z p-value
Fixed effect model   0.1808 [0.0891; 0.2725] 3.87  0.0001
Random effects model 0.2909 [0.0458; 0.5360] 2.33  0.0200

Quantifying heterogeneity:
tau^2 = 0.0640; H = 2.10 [1.41; 3.13]; I^2 = 77.3% [49.6%; 89.8%]

Test of heterogeneity:
     Q d.f. p-value
 22.05    5  0.0005

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=rm.nebf$studlab,pval=rm.nebf$pval))

No USA(CA_MA_MO) data

rmshare.sum.6.nohav<-rmshare.sum.6[!rownames(rmshare.sum.6) %in% "hav",]
rm.nebf<-metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=rmshare.sum.6.nohav,sm="RD", backtransf=FALSE)
forest(rm.nebf,smlab="Standardized \n microbiome age difference",sortvar=rmshare.sum.6.nohav$pop)

rm.nebf

                                           RD            95%-CI %W(fixed)
Thompson et al 2015 ( USA(NC) )        0.8754 [ 0.1366; 1.6142]       1.6
Wood et al 2017 ( South Africa )       0.1214 [-0.2963; 0.5390]       5.0
Pannaraj et al 2017 ( USA(CA_FL) )     0.2665 [ 0.0464; 0.4866]      18.2
Bender et al 2016 ( Haiti )           -0.1769 [-0.8203; 0.4664]       2.1
Subramanian et al 2014 ( Bangladesh )  0.0568 [-0.0605; 0.1741]      63.9
Azad et al 2015 ( Canada )             0.6149 [ 0.3048; 0.9250]       9.1
                                      %W(random)
Thompson et al 2015 ( USA(NC) )              7.4
Wood et al 2017 ( South Africa )            15.0
Pannaraj et al 2017 ( USA(CA_FL) )          22.8
Bender et al 2016 ( Haiti )                  9.0
Subramanian et al 2014 ( Bangladesh )       26.7
Azad et al 2015 ( Canada )                  19.0

Number of studies combined: k = 6

                         RD           95%-CI    z p-value
Fixed effect model   0.1574 [0.0636; 0.2512] 3.29  0.0010
Random effects model 0.2602 [0.0267; 0.4936] 2.18  0.0290

Quantifying heterogeneity:
tau^2 = 0.0495; H = 1.83 [1.20; 2.80]; I^2 = 70.3% [30.5%; 87.3%]

Test of heterogeneity:
     Q d.f. p-value
 16.83    5  0.0048

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=rm.nebf$studlab,pval=rm.nebf$pval))

Vaginal

rm.nebf<-metagen(estimate.nebf, se.nebf, studlab=study,data=rmshare.vag,sm="RD", backtransf=FALSE)
forest(rm.nebf,smlab="Microbiome age difference")

rm.nebf

                                         RD            95%-CI %W(fixed)
Azad et al 2015 (Canada)             2.0223 [ 1.0746; 2.9701]      23.9
Bender et al 2016 (Haiti)           -0.5119 [-2.7851; 1.7614]       4.1
Pannaraj et al 2017 (USA(CA_FL))     0.5934 [-0.0958; 1.2826]      45.1
Sordillo et al 2017 (USA(CA_MA_MO))  1.5749 [ 0.6821; 2.4677]      26.9
                                    %W(random)
Azad et al 2015 (Canada)                  27.7
Bender et al 2016 (Haiti)                 10.8
Pannaraj et al 2017 (USA(CA_FL))          32.8
Sordillo et al 2017 (USA(CA_MA_MO))       28.8

Number of studies combined: k = 4

                         RD           95%-CI    z  p-value
Fixed effect model   1.1523 [0.6894; 1.6152] 4.88 < 0.0001
Random effects model 1.1519 [0.2829; 2.0209] 2.60   0.0094

Quantifying heterogeneity:
tau^2 = 0.4763; H = 1.70 [1.00; 2.92]; I^2 = 65.4% [0.0%; 88.2%]

Test of heterogeneity:
    Q d.f. p-value
 8.68    3  0.0338

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=rm.nebf$studlab,pval=rm.nebf$pval))

C-section

rm.nebf<-metagen(estimate.nebf, se.nebf, studlab=study,data=rmshare.cs,sm="RD", backtransf=FALSE)
forest(rm.nebf,smlab="Microbiome age difference")

rm.nebf

                                         RD            95%-CI %W(fixed)
Azad et al 2015 (Canada)             0.6820 [-1.3391; 2.7031]      14.4
Bender et al 2016 (Haiti)           -0.3508 [-4.3974; 3.6958]       3.6
Pannaraj et al 2017 (USA(CA_FL))     0.7880 [-0.2996; 1.8755]      49.7
Sordillo et al 2017 (USA(CA_MA_MO))  1.1977 [-0.1522; 2.5476]      32.3
                                    %W(random)
Azad et al 2015 (Canada)                  14.4
Bender et al 2016 (Haiti)                  3.6
Pannaraj et al 2017 (USA(CA_FL))          49.7
Sordillo et al 2017 (USA(CA_MA_MO))       32.3

Number of studies combined: k = 4

                         RD           95%-CI    z p-value
Fixed effect model   0.8641 [0.0971; 1.6310] 2.21  0.0272
Random effects model 0.8641 [0.0971; 1.6310] 2.21  0.0272

Quantifying heterogeneity:
tau^2 = 0; H = 1.00 [1.00; 1.17]; I^2 = 0.0% [0.0%; 27.2%]

Test of heterogeneity:
    Q d.f. p-value
 0.63    3  0.8893

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=rm.nebf$studlab,pval=rm.nebf$pval))

Duration (month) of exbf

rmdat.rm<-as.data.frame(rmdat.rm)
b2<-gamm(age.predicted~s(age.sample,by=month.exbf2) +month.exbf2,family=gaussian,
         data=rmdat.rm,random=list(personid=~1))
pred <- predict(b2$gam, newdata = rmdat.rm,se.fit=TRUE)
datfit<-cbind(rmdat.rm, fit=pred$fit,ul=(pred$fit+(1.96*pred$se.fit)),ll=(pred$fit-(1.96*pred$se.fit)))
pexbf2<-ggplot()+ geom_point(data = rmdat.rm, aes(x = age.sample, y = age.predicted, group = personid, colour=month.exbf2))+
  geom_line(data = rmdat.rm, aes(x = age.sample, y = age.predicted, group = personid, colour=month.exbf2),size=0.3)+
  geom_line(data = datfit,aes(x = age.sample, y = fit, colour=month.exbf2),size = 1.5)+
  geom_ribbon(data = datfit,aes(x=age.sample, ymax=ul, ymin=ll,group=month.exbf2), alpha=.5)+
  #theme(legend.position = "bottom")+
  theme(legend.title = element_text(colour="black", size=10))+
  labs(color='Duration EBF')+
  scale_x_continuous(breaks=seq(from=0,to=24,by=3),
                     labels=seq(from=0,to=24,by=3))+
  #theme(legend.text = element_text(colour="black", size = 10))+
  theme(legend.position = c(0.15,0.85),legend.title=element_text(size=8),legend.text=element_text(size=6))+
  theme(legend.key.size = unit(0.5, "cm"),
        axis.line = element_line(colour = "black"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.background = element_blank())+
  xlab("Chronological age (month)") +ylab("Microbiome age (month)")
pexbf2

rmdat.rm.6plus15<-rmdat.rm[rmdat.rm$age.sample>6 & rmdat.rm$age.sample<15,]
b2<-gamm(age.predicted~s(age.sample,by=month.exbf2) +month.exbf2,family=gaussian,
         data=rmdat.rm.6plus15,random=list(personid=~1))

Test for age > 6 months and <15months GAM part

kable(summary(b2$gam)$p.table)

kable(summary(b2$gam)$s.table)

LME part

kable(summary(b2$lme)$tTable)

Chao1

chao1.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="chao1"),sm="RD", backtransf=FALSE)
forest(chao1.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="chao1")$pop)

chao1.nebf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0677 [-0.0136; 0.1491]      79.0
Azad et al 2015 ( Canada )            0.3752 [ 0.0596; 0.6908]       5.3
Bender et al 2016 ( Haiti )           0.4922 [-0.3266; 1.3110]       0.8
Wood et al 2017 ( South Africa )      0.5113 [ 0.0923; 0.9303]       3.0
Pannaraj et al 2017 ( USA(CA_FL) )    0.1201 [-0.1362; 0.3764]       8.0
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.6457 [ 0.2430; 1.0485]       3.2
Thompson et al 2015 ( USA(NC) )       0.4111 [-0.4054; 1.2277]       0.8
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       28.0
Azad et al 2015 ( Canada )                  16.8
Bender et al 2016 ( Haiti )                  4.9
Wood et al 2017 ( South Africa )            12.6
Pannaraj et al 2017 ( USA(CA_FL) )          19.6
Sordillo et al 2017 ( USA(CA_MA_MO) )       13.2
Thompson et al 2015 ( USA(NC) )              4.9

Number of studies combined: k = 7

                         RD           95%-CI    z p-value
Fixed effect model   0.1259 [0.0536; 0.1982] 3.41  0.0006
Random effects model 0.2995 [0.1020; 0.4969] 2.97  0.0030

Quantifying heterogeneity:
tau^2 = 0.0346; H = 1.59 [1.05; 2.41]; I^2 = 60.7% [9.9%; 82.8%]

Test of heterogeneity:
     Q d.f. p-value
 15.25    6  0.0184

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=chao1.nebf$studlab,pval=chao1.nebf$pval))

Observed_species

os.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="observed_species"),sm="RD", backtransf=FALSE)
forest(os.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="observed_species")$pop)

os.nebf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0673 [-0.0127; 0.1473]      77.1
Azad et al 2015 ( Canada )            0.3721 [ 0.0590; 0.6852]       5.0
Bender et al 2016 ( Haiti )           0.4601 [-0.3413; 1.2616]       0.8
Wood et al 2017 ( South Africa )      0.4996 [ 0.0848; 0.9143]       2.9
Pannaraj et al 2017 ( USA(CA_FL) )    0.1447 [-0.0763; 0.3657]      10.1
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.6461 [ 0.2454; 1.0467]       3.1
Thompson et al 2015 ( USA(NC) )       0.3555 [-0.3448; 1.0557]       1.0
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       28.2
Azad et al 2015 ( Canada )                  16.1
Bender et al 2016 ( Haiti )                  4.5
Wood et al 2017 ( South Africa )            12.0
Pannaraj et al 2017 ( USA(CA_FL) )          21.0
Sordillo et al 2017 ( USA(CA_MA_MO) )       12.5
Thompson et al 2015 ( USA(NC) )              5.7

Number of studies combined: k = 7

                         RD           95%-CI    z p-value
Fixed effect model   0.1266 [0.0564; 0.1969] 3.53  0.0004
Random effects model 0.2909 [0.1053; 0.4765] 3.07  0.0021

Quantifying heterogeneity:
tau^2 = 0.0301; H = 1.59 [1.05; 2.41]; I^2 = 60.4% [9.1%; 82.7%]

Test of heterogeneity:
     Q d.f. p-value
 15.14    6  0.0192

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=os.nebf$studlab,pval=os.nebf$pval))

Pd_whole_tree

pwt.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="pd_whole_tree"),sm="RD", backtransf=FALSE)
forest(pwt.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="pd_whole_tree")$pop)

pwt.nebf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0953 [ 0.0196; 0.1710]      80.5
Azad et al 2015 ( Canada )            0.2877 [-0.0327; 0.6082]       4.5
Bender et al 2016 ( Haiti )           0.3082 [-0.3766; 0.9931]       1.0
Wood et al 2017 ( South Africa )      0.4724 [ 0.0489; 0.8958]       2.6
Pannaraj et al 2017 ( USA(CA_FL) )    0.2168 [ 0.0121; 0.4215]      11.0
Thompson et al 2015 ( USA(NC) )       0.7261 [-0.2385; 1.6906]       0.5
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       52.8
Azad et al 2015 ( Canada )                  12.0
Bender et al 2016 ( Haiti )                  3.1
Wood et al 2017 ( South Africa )             7.4
Pannaraj et al 2017 ( USA(CA_FL) )          23.2
Thompson et al 2015 ( USA(NC) )              1.6

Number of studies combined: k = 6

                         RD           95%-CI    z p-value
Fixed effect model   0.1322 [0.0643; 0.2001] 3.82  0.0001
Random effects model 0.1910 [0.0685; 0.3134] 3.06  0.0022

Quantifying heterogeneity:
tau^2 = 0.0059; H = 1.15 [1.00; 1.78]; I^2 = 25.0% [0.0%; 68.4%]

Test of heterogeneity:
    Q d.f. p-value
 6.66    5  0.2469

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=pwt.nebf$studlab,pval=pwt.nebf$pval))

Shannon

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="shannon"),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="shannon")$pop)

shannon.nebf

                                           RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh )  0.2592 [ 0.1185; 0.3999]      57.3
Azad et al 2015 ( Canada )             0.3262 [ 0.0159; 0.6365]      11.8
Bender et al 2016 ( Haiti )           -0.1146 [-0.7954; 0.5662]       2.4
Wood et al 2017 ( South Africa )       0.3071 [-0.1311; 0.7452]       5.9
Pannaraj et al 2017 ( USA(CA_FL) )     0.3732 [ 0.0808; 0.6657]      13.3
Sordillo et al 2017 ( USA(CA_MA_MO) )  0.7684 [ 0.3821; 1.1546]       7.6
Thompson et al 2015 ( USA(NC) )        0.3001 [-0.5308; 1.1310]       1.6
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       40.5
Azad et al 2015 ( Canada )                  15.7
Bender et al 2016 ( Haiti )                  4.0
Wood et al 2017 ( South Africa )             8.9
Pannaraj et al 2017 ( USA(CA_FL) )          17.2
Sordillo et al 2017 ( USA(CA_MA_MO) )       11.0
Thompson et al 2015 ( USA(NC) )              2.7

Number of studies combined: k = 7

                         RD           95%-CI    z  p-value
Fixed effect model   0.3153 [0.2088; 0.4219] 5.80 < 0.0001
Random effects model 0.3359 [0.1959; 0.4759] 4.70 < 0.0001

Quantifying heterogeneity:
tau^2 = 0.0074; H = 1.12 [1.00; 1.67]; I^2 = 20.9% [0.0%; 64.3%]

Test of heterogeneity:
    Q d.f. p-value
 7.59    6  0.2700

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

No Haiti data

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="shannon"&(pop!="Haiti")),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="shannon"&(pop!="Haiti"))$pop)

shannon.nebf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.2592 [ 0.1185; 0.3999]      58.8
Azad et al 2015 ( Canada )            0.3262 [ 0.0159; 0.6365]      12.1
Wood et al 2017 ( South Africa )      0.3071 [-0.1311; 0.7452]       6.1
Pannaraj et al 2017 ( USA(CA_FL) )    0.3732 [ 0.0808; 0.6657]      13.6
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.7684 [ 0.3821; 1.1546]       7.8
Thompson et al 2015 ( USA(NC) )       0.3001 [-0.5308; 1.1310]       1.7
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       45.8
Azad et al 2015 ( Canada )                  15.5
Wood et al 2017 ( South Africa )             8.5
Pannaraj et al 2017 ( USA(CA_FL) )          17.1
Sordillo et al 2017 ( USA(CA_MA_MO) )       10.6
Thompson et al 2015 ( USA(NC) )              2.5

Number of studies combined: k = 6

                         RD           95%-CI    z  p-value
Fixed effect model   0.3261 [0.2183; 0.4340] 5.93 < 0.0001
Random effects model 0.3483 [0.2145; 0.4821] 5.10 < 0.0001

Quantifying heterogeneity:
tau^2 = 0.0050; H = 1.10 [1.00; 1.61]; I^2 = 16.9% [0.0%; 61.6%]

Test of heterogeneity:
    Q d.f. p-value
 6.02    5  0.3047

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

No UNC data

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="shannon"&(pop!="USA(NC)")),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="shannon"&(pop!="USA(NC)"))$pop)

shannon.nebf

                                           RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh )  0.2592 [ 0.1185; 0.3999]      58.3
Azad et al 2015 ( Canada )             0.3262 [ 0.0159; 0.6365]      12.0
Bender et al 2016 ( Haiti )           -0.1146 [-0.7954; 0.5662]       2.5
Wood et al 2017 ( South Africa )       0.3071 [-0.1311; 0.7452]       6.0
Pannaraj et al 2017 ( USA(CA_FL) )     0.3732 [ 0.0808; 0.6657]      13.5
Sordillo et al 2017 ( USA(CA_MA_MO) )  0.7684 [ 0.3821; 1.1546]       7.7
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       36.5
Azad et al 2015 ( Canada )                  17.2
Bender et al 2016 ( Haiti )                  4.9
Wood et al 2017 ( South Africa )            10.3
Pannaraj et al 2017 ( USA(CA_FL) )          18.6
Sordillo et al 2017 ( USA(CA_MA_MO) )       12.6

Number of studies combined: k = 6

                         RD           95%-CI    z  p-value
Fixed effect model   0.3156 [0.2082; 0.4230] 5.76 < 0.0001
Random effects model 0.3427 [0.1851; 0.5004] 4.26 < 0.0001

Quantifying heterogeneity:
tau^2 = 0.0126; H = 1.23 [1.00; 1.94]; I^2 = 34.1% [0.0%; 73.5%]

Test of heterogeneity:
    Q d.f. p-value
 7.58    5  0.1807

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

No USA(CA_MA_MO) data

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="shannon"&(pop!="USA(CA_MA_MO)")),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="shannon"&(pop!="USA(CA_MA_MO)"))$pop)

shannon.nebf

                                           RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh )  0.2592 [ 0.1185; 0.3999]      62.0
Azad et al 2015 ( Canada )             0.3262 [ 0.0159; 0.6365]      12.8
Bender et al 2016 ( Haiti )           -0.1146 [-0.7954; 0.5662]       2.7
Wood et al 2017 ( South Africa )       0.3071 [-0.1311; 0.7452]       6.4
Pannaraj et al 2017 ( USA(CA_FL) )     0.3732 [ 0.0808; 0.6657]      14.4
Thompson et al 2015 ( USA(NC) )        0.3001 [-0.5308; 1.1310]       1.8
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       62.0
Azad et al 2015 ( Canada )                  12.8
Bender et al 2016 ( Haiti )                  2.7
Wood et al 2017 ( South Africa )             6.4
Pannaraj et al 2017 ( USA(CA_FL) )          14.4
Thompson et al 2015 ( USA(NC) )              1.8

Number of studies combined: k = 6

                         RD           95%-CI    z  p-value
Fixed effect model   0.2780 [0.1672; 0.3889] 4.92 < 0.0001
Random effects model 0.2780 [0.1672; 0.3889] 4.92 < 0.0001

Quantifying heterogeneity:
tau^2 = 0; H = 1.00 [1.00; 1.21]; I^2 = 0.0% [0.0%; 32.0%]

Test of heterogeneity:
    Q d.f. p-value
 1.87    5  0.8674

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

Vaginal

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=study,data=subset(alphasum.vag, index=="shannon"),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference")

shannon.nebf

                                         RD            95%-CI %W(fixed)
 Azad et al 2015 (Canada)            0.1424 [-0.0984; 0.3832]      40.8
Bender et al 2016 (Haiti)           -0.2018 [-0.9752; 0.5715]       4.0
Pannaraj et al 2017 (USA(CA_FL))     0.3734 [ 0.1294; 0.6174]      39.7
Sordillo et al 2017 (USA(CA_MA_MO))  0.4321 [ 0.0420; 0.8223]      15.5
                                    %W(random)
 Azad et al 2015 (Canada)                 38.3
Bender et al 2016 (Haiti)                  5.4
Pannaraj et al 2017 (USA(CA_FL))          37.6
Sordillo et al 2017 (USA(CA_MA_MO))       18.6

Number of studies combined: k = 4

                         RD           95%-CI    z p-value
Fixed effect model   0.2656 [0.1118; 0.4194] 3.38  0.0007
Random effects model 0.2647 [0.0798; 0.4495] 2.81  0.0050

Quantifying heterogeneity:
tau^2 = 0.0081; H = 1.13 [1.00; 2.90]; I^2 = 22.3% [0.0%; 88.1%]

Test of heterogeneity:
    Q d.f. p-value
 3.86    3  0.2771

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

C-section

shannon.nebf <- metagen(estimate.nebf, se.nebf, studlab=study,data=subset(alphasum.cs, index=="shannon"),sm="RD", backtransf=FALSE)
forest(shannon.nebf,smlab="Standardized \n diversity difference")

shannon.nebf

                                         RD            95%-CI %W(fixed)
 Azad et al 2015 (Canada)            0.4873 [ 0.0230; 0.9515]      40.3
Pannaraj et al 2017 (USA(CA_FL))    -0.0737 [-0.5869; 0.4395]      33.0
Sordillo et al 2017 (USA(CA_MA_MO))  1.0833 [ 0.5127; 1.6540]      26.7
                                    %W(random)
 Azad et al 2015 (Canada)                 34.9
Pannaraj et al 2017 (USA(CA_FL))          33.4
Sordillo et al 2017 (USA(CA_MA_MO))       31.7

Number of studies combined: k = 3

                         RD            95%-CI    z p-value
Fixed effect model   0.4612 [ 0.1665; 0.7560] 3.07  0.0022
Random effects model 0.4888 [-0.1327; 1.1103] 1.54  0.1232

Quantifying heterogeneity:
tau^2 = 0.2323; H = 2.09 [1.16; 3.77]; I^2 = 77.1% [25.6%; 93.0%]

Test of heterogeneity:
    Q d.f. p-value
 8.75    2  0.0126

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.nebf$studlab,pval=shannon.nebf$pval))

Put meta-analysis results (random models) of all indexes together

a.nebf<-ggplot(data=a.metatab.r,aes(x=estimate.nebf,y=index))+
  geom_point(shape=17, colour="red")+
  geom_errorbarh(aes(xmin=ll.nebf,xmax=ul.nebf),height=0.0, colour="red")+
  scale_x_continuous(breaks=seq(from=0,to=1,by=0.2),
                     labels=seq(from=0,to=1,by=0.2))+
  geom_vline(xintercept=0,linetype="dashed")+
  theme(axis.line = element_line(colour = "black"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.background = element_blank())+
  ggtitle("Non-EBF vs. EBF")+
  xlab("Pooled Standardized diversity difference")+ylab("Alpha diversity index")
a.nebf

kable(a.metatab.r[,grep(".nebf",colnames(a.metatab.r))])

Chao1

chao1.conbf <- metagen(estimate.conbf, se.conbf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="chao1"),sm="RD", backtransf=FALSE)
forest(chao1.conbf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="chao1")$pop)

chao1.conbf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0550 [-0.0189; 0.1289]      65.2
Azad et al 2015 ( Canada )            0.5725 [ 0.3986; 0.7464]      11.8
Bender et al 2016 ( Haiti )               NA                         0.0
Wood et al 2017 ( South Africa )          NA                         0.0
Pannaraj et al 2017 ( USA(CA_FL) )    0.1302 [-0.0599; 0.3202]       9.9
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.4245 [ 0.2553; 0.5937]      12.4
Thompson et al 2015 ( USA(NC) )       0.4233 [-0.2687; 1.1153]       0.7
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       24.7
Azad et al 2015 ( Canada )                  22.4
Bender et al 2016 ( Haiti )                  0.0
Wood et al 2017 ( South Africa )             0.0
Pannaraj et al 2017 ( USA(CA_FL) )          21.9
Sordillo et al 2017 ( USA(CA_MA_MO) )       22.6
Thompson et al 2015 ( USA(NC) )              8.3

Number of studies combined: k = 5

                         RD           95%-CI    z  p-value
Fixed effect model   0.1720 [0.1123; 0.2316] 5.65 < 0.0001
Random effects model 0.3016 [0.0580; 0.5451] 2.43   0.0152

Quantifying heterogeneity:
tau^2 = 0.0610; H = 3.13 [2.19; 4.49]; I^2 = 89.8% [79.1%; 95.0%]

Test of heterogeneity:
     Q d.f.  p-value
 39.26    4 < 0.0001

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=chao1.conbf$studlab,pval=chao1.conbf$pval))

Observed_species

os.conbf <- metagen(estimate.conbf, se.conbf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="observed_species"),sm="RD", backtransf=FALSE)
forest(os.conbf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="observed_species")$pop)

os.conbf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0603 [-0.0124; 0.1329]      63.5
Azad et al 2015 ( Canada )            0.5895 [ 0.4168; 0.7623]      11.2
Bender et al 2016 ( Haiti )               NA                         0.0
Wood et al 2017 ( South Africa )          NA                         0.0
Pannaraj et al 2017 ( USA(CA_FL) )    0.1163 [-0.0479; 0.2806]      12.4
Sordillo et al 2017 ( USA(CA_MA_MO) ) 0.4610 [ 0.2928; 0.6291]      11.9
Thompson et al 2015 ( USA(NC) )       0.3661 [-0.2275; 0.9598]       1.0
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       24.1
Azad et al 2015 ( Canada )                  21.9
Bender et al 2016 ( Haiti )                  0.0
Wood et al 2017 ( South Africa )             0.0
Pannaraj et al 2017 ( USA(CA_FL) )          22.1
Sordillo et al 2017 ( USA(CA_MA_MO) )       22.0
Thompson et al 2015 ( USA(NC) )             10.0

Number of studies combined: k = 5

                         RD           95%-CI    z  p-value
Fixed effect model   0.1772 [0.1192; 0.2351] 6.00 < 0.0001
Random effects model 0.3070 [0.0641; 0.5500] 2.48   0.0133

Quantifying heterogeneity:
tau^2 = 0.0625; H = 3.30 [2.33; 4.69]; I^2 = 90.8% [81.6%; 95.4%]

Test of heterogeneity:
     Q d.f.  p-value
 43.68    4 < 0.0001

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=os.conbf$studlab,pval=os.conbf$pval))

Pd_whole_tree

pwt.conbf <- metagen(estimate.conbf, se.conbf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="pd_whole_tree"),sm="RD", backtransf=FALSE)
forest(pwt.conbf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="pd_whole_tree")$pop)

pwt.conbf

                                          RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh ) 0.0646 [-0.0047; 0.1338]      72.9
Azad et al 2015 ( Canada )            0.5311 [ 0.3540; 0.7082]      11.2
Bender et al 2016 ( Haiti )               NA                         0.0
Wood et al 2017 ( South Africa )          NA                         0.0
Pannaraj et al 2017 ( USA(CA_FL) )    0.2032 [ 0.0523; 0.3541]      15.4
Thompson et al 2015 ( USA(NC) )       0.8566 [ 0.0509; 1.6623]       0.5
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       32.9
Azad et al 2015 ( Canada )                  29.1
Bender et al 2016 ( Haiti )                  0.0
Wood et al 2017 ( South Africa )             0.0
Pannaraj et al 2017 ( USA(CA_FL) )          30.2
Thompson et al 2015 ( USA(NC) )              7.8

Number of studies combined: k = 4

                         RD           95%-CI    z  p-value
Fixed effect model   0.1422 [0.0830; 0.2013] 4.71 < 0.0001
Random effects model 0.3037 [0.0475; 0.5598] 2.32   0.0202

Quantifying heterogeneity:
tau^2 = 0.0506; H = 3.00 [1.97; 4.57]; I^2 = 88.9% [74.2%; 95.2%]

Test of heterogeneity:
     Q d.f.  p-value
 26.99    3 < 0.0001

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=pwt.conbf$studlab,pval=pwt.conbf$pval))

Shannon

shannon.conbf <- metagen(estimate.conbf, se.conbf, studlab=paste(author,year,"(",pop,")"),data=subset(alphaall,index=="shannon"&(!pop %in% c("Haiti","South Africa"))),sm="RD", backtransf=FALSE)
forest(shannon.conbf,smlab="Standardized \n diversity difference",sortvar=subset(alphaall,index=="shannon"&(!pop %in% c("Haiti","South Africa")))$pop)

shannon.conbf

                                           RD            95%-CI %W(fixed)
Subramanian et al 2014 ( Bangladesh )  0.2131 [ 0.0843; 0.3420]      39.0
Azad et al 2015 ( Canada )             0.5872 [ 0.4153; 0.7590]      21.9
Pannaraj et al 2017 ( USA(CA_FL) )     0.2570 [ 0.0375; 0.4765]      13.5
Sordillo et al 2017 ( USA(CA_MA_MO) )  0.5660 [ 0.4038; 0.7282]      24.6
Thompson et al 2015 ( USA(NC) )       -0.0180 [-0.8408; 0.8049]       1.0
                                      %W(random)
Subramanian et al 2014 ( Bangladesh )       25.8
Azad et al 2015 ( Canada )                  23.8
Pannaraj et al 2017 ( USA(CA_FL) )          21.4
Sordillo et al 2017 ( USA(CA_MA_MO) )       24.2
Thompson et al 2015 ( USA(NC) )              4.8

Number of studies combined: k = 5

                         RD           95%-CI    z  p-value
Fixed effect model   0.3858 [0.3053; 0.4663] 9.39 < 0.0001
Random effects model 0.3857 [0.1872; 0.5843] 3.81   0.0001

Quantifying heterogeneity:
tau^2 = 0.0355; H = 2.19 [1.42; 3.37]; I^2 = 79.1% [50.4%; 91.2%]

Test of heterogeneity:
     Q d.f. p-value
 19.16    4  0.0007

Details on meta-analytical method:
- Inverse variance method
- DerSimonian-Laird estimator for tau^2

kable(cbind(study=shannon.conbf$studlab,pval=shannon.conbf$pval))

Put meta-analysis results (random models) of all indexes together

a.conbf<-ggplot(data=a.metatab.r,aes(x=estimate.conbf,y=index))+
  geom_point(shape=17, colour="red")+
  geom_errorbarh(aes(xmin=ll.conbf,xmax=ul.conbf),height=0.0, colour="red")+
  geom_vline(xintercept=0,linetype="dashed")+
  scale_x_continuous(breaks=seq(from=0,to=1,by=0.2),
                     labels=seq(from=0,to=1,by=0.2))+
  theme(axis.line = element_line(colour = "black"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.background = element_blank())+
  ggtitle("Trend across EBF, non-EBF, non-BF")+
  xlab("Pooled standardized diversity difference")+ylab("Alpha diversity index")
a.conbf

kable(a.metatab.r[,grep(".conbf",colnames(a.metatab.r))])

l2

kable(metatab.show(metatab=metatab.zi.nounc$random,taxacom.pooled.tab=taxacom.zi.nounc,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l4

kable(metatab.show(metatab=metatab.zi.nounc$random,taxacom.pooled.tab=taxacom.zi.nounc,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l5

kable(metatab.show(metatab=metatab.zi.nounc$random,taxacom.pooled.tab=taxacom.zi.nounc,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l6

kable(metatab.show(metatab=metatab.zi.nounc$random,taxacom.pooled.tab=taxacom.zi.nounc,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l2

kable(metatab.show(metatab=metatab.zi.noha$random,taxacom.pooled.tab=taxacom.zi.noha,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l4

kable(metatab.show(metatab=metatab.zi.noha$random,taxacom.pooled.tab=taxacom.zi.noha,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l5

kable(metatab.show(metatab=metatab.zi.noha$random,taxacom.pooled.tab=taxacom.zi.noha,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l6

kable(metatab.show(metatab=metatab.zi.noha$random,taxacom.pooled.tab=taxacom.zi.noha,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l2

kable(metatab.show(metatab=metatab.zi.nohav$random,taxacom.pooled.tab=taxacom.zi.nohav,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l4

kable(metatab.show(metatab=metatab.zi.nohav$random,taxacom.pooled.tab=taxacom.zi.nohav,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l5

kable(metatab.show(metatab=metatab.zi.nohav$random,taxacom.pooled.tab=taxacom.zi.nohav,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l6

kable(metatab.show(metatab=metatab.zi.nohav$random,taxacom.pooled.tab=taxacom.zi.nohav,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

l2

kable(metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="main",phyla.col="rainbow",leg.key.size=0.5,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l4

kable(metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l5

kable(metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l6

kable(metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.vag$random,taxacom.pooled.tab=taxacom.zi.vag,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l2

kable(metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="main",phyla.col="rainbow",leg.key.size=0.4,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l4

kable(metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l4",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l5

kable(metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l5",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

l6

kable(metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=metatab.zi.cs$random,taxacom.pooled.tab=taxacom.zi.cs,tax.lev="l6",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="taxa",level="sub",phyla.col="rainbow",leg.key.size=1,leg.text.size=8,heat.text.x.size=6,forest.axis.text.y=8,forest.axis.text.x=8)

No Haiti data

Level 2

kable(metatab.show(metatab=pathmetatab.zi.noha$random,taxacom.pooled.tab=pathcom.zi.noha,sumvar="path",tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=pathmetatab.zi.noha$random,taxacom.pooled.tab=pathcom.zi.noha,sumvar="path",tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="path",leg.key.size=1,leg.text.size=8,heat.text.x.size=4,forest.axis.text.y=7,forest.axis.text.x=7)

Level 3

kable(metatab.show(metatab=pathmetatab.zi.noha$random,taxacom.pooled.tab=pathcom.zi.noha,sumvar="path",tax.lev="l3",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=pathmetatab.zi.noha$random,taxacom.pooled.tab=pathcom.zi.noha,sumvar="path",tax.lev="l3",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="path",leg.key.size=1,leg.text.size=8,heat.text.x.size=4,forest.axis.text.y=6,forest.axis.text.x=6)

No USA(CA_MA_MO) data

Level 2

kable(metatab.show(metatab=pathmetatab.zi.nohav$random,taxacom.pooled.tab=pathcom.zi.nohav,sumvar="path",tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=pathmetatab.zi.nohav$random,taxacom.pooled.tab=pathcom.zi.nohav,sumvar="path",tax.lev="l2",showvar=".nebf",p.cutoff.type="p", p.cutoff=1,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="path",leg.key.size=1,leg.text.size=8,heat.text.x.size=4,forest.axis.text.y=7,forest.axis.text.x=7)

Level 3

kable(metatab.show(metatab=pathmetatab.zi.nohav$random,taxacom.pooled.tab=pathcom.zi.nohav,sumvar="path",tax.lev="l3",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="table"))

metadat<-metatab.show(metatab=pathmetatab.zi.nohav$random,taxacom.pooled.tab=pathcom.zi.nohav,sumvar="path",tax.lev="l3",showvar=".nebf",p.cutoff.type="p", p.cutoff=0.05,display="data",fill.value="log(OR)")
meta.niceplot(metadat=metadat,sumtype="path",leg.key.size=1,leg.text.size=8,heat.text.x.size=4,forest.axis.text.y=6,forest.axis.text.x=6)

Phylum (L2)

p.exbf.l2<-taxa.mean.plot(tabmean=taxa.meansdn.exbf2.rm,taxlist=taxlist.rm,tax.lev="l2", comvar="month.exbf2", groupvar="age.sample",mean.filter=0.005)
p.exbf.l2$p

#dev.off()

GAMLSS

kable(taxcomtab.show(taxcomtab=taxacom.6plus.exbf2.zi.rm,tax.select=p.exbf.l2$taxuse.rm, tax.lev="l2",p.adjust.method="fdr"))

Order (l4)

p.exbf.l4<-taxa.mean.plot(tabmean=taxa.meansdn.exbf2.rm,taxlist=taxlist.rm,tax.lev="l4", comvar="month.exbf2", groupvar="age.sample",mean.filter=0.005)
p.exbf.l4$p

GAMLSS

kable(taxcomtab.show(taxcomtab=taxacom.6plus.exbf2.zi.rm,tax.select=p.exbf.l4$taxuse.rm, tax.lev="l4",p.adjust.method="fdr"))

Family (L5)

p.exbf.l5<-taxa.mean.plot(tabmean=taxa.meansdn.exbf2.rm,taxlist=taxlist.rm,tax.lev="l5", comvar="month.exbf2", groupvar="age.sample",mean.filter=0.005,show.taxname="short",legend.position="right")
p.exbf.l5$p

# better plot view
for (i in 1: length(taxa.meansdn.exbf2.rm)){
  taxa.meansdn.exbf2.rm[[i]]$month.exbf2l<-mapvalues(taxa.meansdn.exbf2.rm[[i]]$month.exbf2,from=c("<=2 months",">2 months"),to=c("Duration EBF <=2 months","Duration EBF >2 months"))
}
p.exbf.l5.noleg<-taxa.mean.plot(tabmean=taxa.meansdn.exbf2.rm,taxlist=taxlist.rm,tax.lev="l5", comvar="month.exbf2l", groupvar="age.sample",mean.filter=0.005,legend.position="none")
grid.arrange(pexbf2,p.exbf.l5.noleg$p,nrow=1)

GAMLSS

kable(taxcomtab.show(taxcomtab=taxacom.6plus.exbf2.zi.rm,tax.select=p.exbf.l5$taxuse.rm, tax.lev="l5",p.adjust.method="fdr"))

Genus (L6)

p.exbf.l6<-taxa.mean.plot(tabmean=taxa.meansdn.exbf2.rm,taxlist=taxlist.rm,tax.lev="l6", comvar="month.exbf2", groupvar="age.sample",mean.filter=0.005)
p.exbf.l6$p

GAMLSS

kable(taxcomtab.show(taxcomtab=taxacom.6plus.exbf2.zi.rm,tax.select=p.exbf.l6$taxuse.rm, tax.lev="l6",p.adjust.method="fdr"))

Stratified by duration of exclusive bf

p.dia.exbf2.6plus.l5<-taxa.mean.plot(tabmean=taxa.meansdn.dia.exbf2.6plus.rm,taxlist=taxlist.rm,tax.lev="l5", comvar="diarrhea", groupvar="month.exbf2",mean.filter=0.005,legend.position="right")
p.dia.exbf2.6plus.l5$p

#more detail legend
teste<-taxa.meansdn.dia.exbf2.6plus.rm
for (i in 1:length(names(teste))){
  teste[[i]]$month.exbf2l<-mapvalues(teste[[i]]$month.exbf2,from=c("<=2 months",">2 months"),to=c("Duration EBF <=2 months","Duration EBF >2 months"))
}
p.dia.exbf2.6plus.l5noleg<-taxa.mean.plot(tabmean=teste,taxlist=taxlist.rm,tax.lev="l5", comvar="diarrhea", groupvar="month.exbf2l",mean.filter=0.005,legend.position="none",ylab="Relative abundance (6 months - 2 years)")

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.exbf2.zi.rm, tax.lev="l5",tax.select=p.dia.exbf2.6plus.l5$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.exbf2plus.zi.rm, tax.lev="l5",tax.select=p.dia.exbf2.6plus.l5$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

p.dia.exbf2.6plus.l6<-taxa.mean.plot(tabmean=taxa.meansdn.dia.exbf2.6plus.rm,taxlist=taxlist.rm,tax.lev="l6", comvar="diarrhea", groupvar="month.exbf2",mean.filter=0.005,legend.position = "right")
p.dia.exbf2.6plus.l6$p

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.exbf2.zi.rm, tax.lev="l6",tax.select=p.dia.exbf2.6plus.l6$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.exbf2plus.zi.rm, tax.lev="l6",tax.select=p.dia.exbf2.6plus.l6$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

Stratified by bf status

taxa.meansdn.dia.bf.6plus.rm.s<-list()
for (i in 1:5){
  taxa.meansdn.dia.bf.6plus.rm.s[[i]]<-taxa.meansdn.dia.bf.6plus.rm[[i]][taxa.meansdn.dia.bf.6plus.rm[[i]]$bf!="ExclusiveBF",]
  taxa.meansdn.dia.bf.6plus.rm.s[[i]]$bf<-drop.levels(taxa.meansdn.dia.bf.6plus.rm.s[[i]]$bf,reorder=FALSE)
}
names(taxa.meansdn.dia.bf.6plus.rm.s)<-names(taxa.meansdn.dia.bf.6plus.rm)

p.dia.bf.6plus.l5<-taxa.mean.plot(tabmean=taxa.meansdn.dia.bf.6plus.rm.s,tax.select=p.dia.exbf2.6plus.l5$taxuse.rm,taxlist=taxlist.rm,tax.lev="l5", comvar="diarrhea", groupvar="bf",mean.filter=0.005,show.taxname = "short")
p.dia.bf.6plus.l5$p

# reverse levels for better combined plot view
testr=taxa.meansdn.dia.bf.6plus.rm.s
for (i in 1: length(names(testr))){
  testr[[i]]$bfl<-mapvalues(testr[[i]]$bf,from=c("No_BF","Non_exclusiveBF"),to=c("No BF when diarrhea","BF when diarrhea"))
  testr[[i]]$bfl<-factor(testr[[i]]$bfl,levels=c("No BF when diarrhea","BF when diarrhea"))
}
p.dia.bf.6plus.l5.nolegend<-taxa.mean.plot(tabmean=testr,tax.select=p.dia.exbf2.6plus.l5$taxuse.rm,taxlist=taxlist.rm,tax.lev="l5", comvar="diarrhea", groupvar="bfl",mean.filter=0.005,show.taxname = "short",legend.position = "none",ylab="Relative abundance (6 months - 2 years)")

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.nexbf.zi.rm, tax.lev="l5",tax.select=p.dia.bf.6plus.l5$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.nobf.zi.rm, tax.lev="l5",tax.select=p.dia.bf.6plus.l5$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

p.dia.bf.6plus.l6<-taxa.mean.plot(tabmean=taxa.meansdn.dia.bf.6plus.rm.s,tax.select=p.dia.exbf2.6plus.l6$taxuse.rm,taxlist=taxlist.rm,tax.lev="l6", comvar="diarrhea", groupvar="bf",mean.filter=0.005)
p.dia.bf.6plus.l6$p

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.nexbf.zi.rm, tax.lev="l6",tax.select=p.dia.bf.6plus.l6$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))

kable(taxcomtab.show(taxcomtab=taxacom.6plus.dia.nobf.zi.rm, tax.lev="l6",tax.select=p.dia.bf.6plus.l6$taxuse.rm,p.adjust.method="fdr",p.cutoff=0.1))